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1 BASICS 



1 Basics 



(AB)- 1 = B X A 1 (1) 

(ABC...) -1 = ...C^B^A 1 (2) 

(A^ 1 = (A~Y (3) 

(A + B) T = A T + B T (4) 

(AB) T = B T A T (5) 

(ABC...) T = ...C T B T A T (6) 

(A^)- 1 = (A- 1 )" (7) 

(A + B) H = A H + B H (8) 

(AB) H = B H A H (9) 

(ABC...) ff = ...C H B H A H (10) 

1.1 Trace and Determinants 

Tr(A) = (11) 

Tr(A) = £.Ai, A 4 = eig(A) (12) 

Tr(A) = Tr(A T ) (13) 

Tr(AB) = Tr(BA) (14) 

Tr(A + B) = Tr(A) + Tr(B) (15) 

Tr(ABC) = Tr(BCA) = Tr(CAB) (16) 

det(A) = I1A A 4 = eig(A) (17) 

det(cA) = c"det(A), if A e R nxn (18) 

det(A T ) = det(A) (19) 

det(AB) = det(A) det(B) (20) 

det(A" 1 ) = l/det(A) (21) 

det(A") = det(A)" (22) 

det(I + uv T ) = l + u T v (23) 

det(I + eA) ^ l + eTr(A), e small (24) 



1.2 The Special Case 2x2 

Consider the matrix A 

A r a 1x a 12 ' 

I A 2 \ A22 

Determinant and trace 

det(A) = AnA 22 - A 12 A 21 (25) 
Tr(A) = A n + A 22 (26) 
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1.2 The Special Case 2x2 



1 BASICS 



Eigenvalues 



A 2 - A • TV (A) + dct(A) = 



Tr(A) + a/Ti^A) 2 — 4dct(A) Tr(A) - ^/Ti^A) 2 — 4dct(A) 
- 2 A2 - 2 



Ai + A 2 = Tr(A) AiA 2 = dct(A) 



Eigenvectors 



Inverse 





Al2 




A 12 


Vl oc 


Xi-Au _ 


v 2 cx 


A 2 - An 



A- 1 



dct(A) 



A 22 — A12 
-Mi An 



(27) 



Petersen & Pedersen, The Matrix Cookbook, Version: November 14, 2008, Page 6 



2 DERIVATIVES 



2 Derivatives 



This section is covering differentiation of a number of expressions with respect to 
a matrix X. Note that it is always assumed that X has no special structure, i.e. 
that the elements of X are independent (e.g. not symmetric, Toeplitz, positive 
definite). See section 2.8 for differentiation of structured matrices. The basic 



assumptions can be written in a formula as 



dX k 



ki 



9Xi 



SikS, 



ikOlj 



(28) 



that is for e.g. vector forms, 



dy 



dxj 
dy 



dx 
dy 



dx 
dy t 



dy 



dxj 
9yj 



The following rules are general and very useful when deriving the differential of 
an expression ([!§]): 



dA 

d{aX) 
d(X + Y) 
9(Tr(X)) 
d(XY) 
0(X o Y) 
9(X ® Y) 

^(x- 1 ) 

<9(det(X)) 
<9(ln(det(X))) 
dX T 

dx H 





adX 
dX + dY 
Tv(dX) 

(dX)Y + X(dY) 

(dX) o Y + X o (<9Y) 

(dX) ® Y + X ® (9Y) 

-X^dXjX" 1 

det(X)Tr(X" 1 9X) 

T^X^SX) 

(dXf 

(dX) H 



(A is a constant) 



(29) 
(30) 
(31) 
(32) 
(33) 
(34) 
(35) 
(36) 
(37) 
(38) 
(39) 
(40) 



2.1 Derivatives of a Determinant 
2.1.1 General form 



ddet(Y) 

dx 

a ddct(Y) 
° dx 

dx 



= dct(Y)Tr 



r-l3Y 

dx 



det(Y) 



+Tr 



-Tr 



Tr 



r — 1 dx 



dx 



_i<9Y 

dx 



Tr 



(41) 



_!<9Y 

dx 



(42) 
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2.2 Derivatives of an Inverse 



2 DERIVATIVES 



2.1.2 Linear forms 

5d6t(X) = det(X),X V ,J3, 



ax 

E^ 1 ^ = ^' det W ( 44 ) 

k lk 

9det(AXB) = det(AXB)(x -i ) T = dct(AXB)(x T ) -i (45) 

2.1.3 Square forms 

If X is square and invertible, then 
<9det(X T AX 



OX 

If X is not square but A is symmetric, then 



2det(X T AX)X~ T (46) 



^ d6t( 5 AX) = 2det(X T AX)AX(X T AX)- 1 



(47) 



If X is not square and A is not symmetric, then 

ddet(XjAX) = dct(x T AX)(AX(x r AX) -i + A T x(x T A T xr i } (48) 
aX 

2.1.4 Other nonlinear forms 

Some special cases are (See [HI [7]) 

<91ndet(X T X)| 



ax 

<91ndet(X T X) _ oyT 
<9X+ 
51n|det(X)| 



2(X+) T (49) 
2X T (50) 



ax 

9det(X fc ) 



ax 

2.2 Derivatives of an Inverse 

From [27J we have the basic identity 

ay- 1 _ v _!aY 



= (X- 1 f = (X r )- 1 (51) 
= fcdet(X fc )X- T (52) 



—Y- 1 !^ 1 (53) 



dx dx 
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2.3 Derivatives of Eigenvalues 



2 DERIVATIVES 



from which it follows 



dX i3 
9a r X- 1 b 

ax 

9det(X- 1 ) 
dX 

<9Tr(AX _1 B) 
dX 

aTr((X + A)- 1 ) 
<9X 



-pt-^wcx- 1 ),-, 

-X- T ab T X- T 
-det(X- 1 )(X- 1 ) T 

-(x-'bax-y 

-((X + A)- 1 (X + A)- 1 ) T 



(54) 

(55) 
(56) 
(57) 
(58) 



From [32 wc have the following result: Let A be an n x n invcrtible square 
matrix, W be the inverse of A, and J(A) is an n x n -variate and diffcrcntiable 
function with respect to A, then the partial differentials of J with respect to A 
and W satisfy 

dJ . _ T d,J . _rp 

— A ' A ' 

dA dW 



2.3 Derivatives of Eigenvalues 

E ci 8'( X ) = 7^ Tr ( X ) = 1 



dX 

d 



J]eig(X) 



ax 

8 



det(X) 



det(X)X 



— T 



2.4 Derivatives of Matrices, Vectors and Scalar Forms 
2.4.1 First Order 



(59) 
(60) 



5x T a 
<9x 
c>a T Xb 
dX 
8a T X T b 
dX 
da T Xa 

dX 
OX 

9(XA) i3 - 
d(X T A) tJ 



da T y 
<9x 

ab T 



ba T 

da T X T a 

ax 



= a 



= aa 



(J mn A)« 
(J nm A) l: , 



(61) 
(62) 
(63) 
(64) 
(65) 

(66) 
(67) 
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2.4 Derivatives of Matrices, Vectors and Scalar Forms 2 DERIVATIVES 



2.4.2 Second Order 

XkiX mn — 2^JXy (68) 

zlmn kl 

db T X T Xc 
dX 

.(Bx + bfCCDx + d) = B T c(Dx + d) + D T G T (Bx + b) (70) 



ax 

J klmn 

fluTvTv„ 

X(bc T + cb T ) (69) 



<9x 
<9(X T BX 



0X« 



S l3 (X T B) kl + S k3 (BX) u (71) 



flfX^RX"! 

' = X'B.J J BX (J ij ) k i = 6^ (72) 



dX l3 



See Sec |9.7| for useful properties of the Single-entry matrix J lJ 

(B + B T )x (73) 



<9x T Bx 
<9x 

ab T x T DXc 



D T Xbc T + DXcb T (74) 



9X 



^(Xb + c) T D(Xb + c) = (D + D T )(Xb + c)b T (75) 
(/X 

Assume W is symmetric, then 

^-(x- As) T W(x- As) = -2A T W(x-As) (76) 
as 

-|-(x-s) T W(x-s) = 2W(x-s) (77) 

J^(x-s) T W(x-s) = -2W(x-s) (78) 

J^(x- As) T W(x- As) = 2W(x-As) (79) 

— (x- As) T W(x- As) = -2W(x-As)s T (80) 
dA 

As a case with complex values the following holds 

d (a ~^ b)2 = -2b(a-x"b)* (81) 
This formula is also known from the LMS algorithm [2] 
2.4.3 Higher order and non-linear 



<9(X n ) 



ra-l 



A7 



J r—0 



^( X rjii X n-l-rj w (82) 
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2.5 Derivatives of Traces 



2 DERIVATIVES 



For proof of the above, see B.l.l 
d 



ax 



a T X"b = ^(X'') T ab T (X"- 1 -'') T 



(83) 



ax 



n—x 

a T (X") T X"b = [x n - 1 - r ab T (X") T X r 



r=0 



.(X r )' i 'X n ab T (X n - 1 - r ) T 



(84) 



See |B.1.1| for a proof. 

Assume s and r are functions of x, i.e. s = s(x),r = r(x), and that A is a 
constant, then 



d_ 
<9x 



s T Ar 



~ds~ 


T 


dr 




Ar + 









A J s 



(85) 



d (Ax) i (Ax) d x J A J Ax 



dx (Bx) T (Bx) ~ 9x x T B T Bx 

A T Ax x T A T AxB T Bx 
x T BBx (x T B T Bx) 2 

2.4.4 Gradient and Hessian 

Using the above we have for the gradient and the Hessian 

/ = x T Ax + b T x 



df 
9x 



d 2 f 



(A + A T )x + b 



= A + A 



(86) 
(87) 



(88) 
(89) 

(90) 



dxdx T 

2.5 Derivatives of Traces 

Assume F(K) to be a diffcrcntiablc function of each of the elements of X. It 
then holds that 

dTr(F(X)) _ T 

dx - /(X) 

where /(•) is the scalar derivative of F(-). 
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2.5 Derivatives of Traces 



2 DERIVATIVES 



2.5.1 First Order 



<9X 



Tr(X) 



J^(AXB) 
9 Tr(AX T B) 



dX 



Tr(X T A) 



d_ 
dX 

3 Tr(A®X) 



ax 



2.5.2 Second Order 



dX 



dX 



d_ 
dX 
_9 
dX 

d 



dX 



d_ 
dX 



dX 



dX 



dX 



Tr 



ax 



I 

A T 

A T B T 
BA 
A 
A 

Tr(A)I 



Tr(X 2 ) = 



Tr(X 2 B) = 



Tr(X T BX) = 
Tr(XBX T ) = 



Tr(AXBX) = 



Tr(X T X) = 



Tr(BXX T ) = 



Tr(B T X T CXB) = 
A T r[X T BXC] = 
^Tr(AXBX T C) = 

(AXB + C)(AXC + C) T 

d 



2X T 

(XB + BX) T 
BX + B T X 
XB T + XB 

A T X T B T + B T X T A T 

2X 

(B + B T )X 
C T XBB T + CXBB T 
BXC + B T XC T 
A T C T XB T + CAXB 



Tr(X <8> X) = 



2A T (AXB + C)B T 

_d_ 
dX 



(91) 
(92) 
(93) 
(94) 
(95) 
(96) 
(97) 



(98) 
(99) 
(100) 
(101) 
(102) 
(103) 
(104) 
(105) 
(106) 
(107) 
(108) 

Tr(X)Tr(X) = 2Tr(X)I(109) 



See 0. 
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2.6 Derivatives of vector norms 



2 DERIVATIVES 



2.5.3 Higher Order 



ax 



Tr(X fc ) 



■k-l\T 



OX 



Tr(AX l 



fc(X fc - i ) 

fc-i 
^(X r AX 



fc— r-l\T 



(110) 

(111) 



r=0 

^Tr [B T X T CXX T CXB] = CXX T CXBB T 

+C T XBB T X T C T X 

+CXBB T X T CX 

+C T XX T C T XBB T 



2.5.4 Other 

d 



dX 



^(AX"^) = -(X^BAX-y = -X~ T A T B T X 



(112) 



(113) 



Assume B and C to be symmetric, then 

a 



dX 



dX 



d_ 
dX 



Tr 



Tr 



Tr (X T CX)" 1 A = -(CX(X T CX)" 1 )(A + A T )(X T CX)" 1 (114) 

(X T CX)- 1 (X T BX)j = -2CX(X T CX)- 1 X T BX(X T CX)- 1 

+2BX(X T CX)" 1 (115) 
(A + X T CX)- 1 (X T BX)j = -2CX(A + X T CX)- 1 X T BX(A + X T CX)- 1 

+2BX(A + X T CX)" 1 (116) 



See [7]. 



«9Tr(sin(X)) 
dX 



= cos(X) rj 



(117) 



2.6 Derivatives of vector norms 
2.6.1 Two-norm 



d_ 

d x a 

<9x ||x - a|| 2 



x — a 



|x-a|| 2 
I 



l x - a \\2 

(x — a)(x — a) T 



l x - a \\2 



|x T x|| 2 



<9x 9x 

2.7 Derivatives of matrix norms 

For more on matrix norms, see Sec. |104| 



x — a 



2x 



(118) 
(119) 
(120) 
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2.8 Derivatives of Structured Matrices 



2 DERIVATIVES 



2.7.1 Probenius norm 

d 



,,v II 2 



dX 



Tr(XX ff ) = 2X 



(121) 



See (2271. Note that this is also a special case of the result in equation 108 



2.8 Derivatives of Structured Matrices 

Assume that the matrix A has some structure, i.e. symmetric, toeplitz, etc. 
In that case the derivatives of the previous section does not apply in general. 
Instead, consider the following general rule for differentiating a scalar function 
/(A) 



df OA 



H 9Akl 9A 



= Tr 



dA 



OA 

dA~< 



(122) 



The matrix differentiated with respect to itself is in this document referred to 
as the structure matrix of A and is defined simply by 



dA 



= S l 



(123) 



If A has no special structure we have simply S lJ = J 1 -' , that is, the structure 
matrix is simply the singleentry matrix. Many structures have a representation 
in singleentry matrices, see Sec. |9.7.6| for more examples of structure matrices. 

2.8.1 The Chain Rule 

Sometimes the objective is to find the derivative of a matrix which is a function 
of another matrix. Let U = /(X), the goal is to find the derivative of the 
function g(U) with respect to X: 



dg(XJ) = 0g(/(X)) 
dX dX 

Then the Chain Rule can then be written the following way: 

M N 



dg(TJ) _ dg(U) 



dX 



dxi 



EE 



k=l 1=1 

Using matrix notation, this can be written as: 



dg(U) 

dXu 



= Tr 



dU 1 OX, 



dg(U) du H 

dun dxij 



d\5 



(124) 



(125) 



(126) 
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2.8 Derivatives of Structured Matrices 



2 DERIVATIVES 



2.8.2 Symmetric 

If A is symmetric, then S lJ 



dA 



That is, e.g., (0): 



: jy + ;p' - J .1 and therefore 



<9Tr(AX) 

<9X 
9det(X) 
<9X 
dlndet(X) 

ax 

2.8.3 Diagonal 

If X is diagonal, then (|19j): 





+ 




T 

— diag 






8A_ 




dA_ 


9A 






A + A T 


-(Aol) 


see 


(131 




det(X)(2X _1 - (X" 1 o 


I)) 




2X 


-l 


(X^ol) 







5Tr(AX) 

ax 



= Aol 



2.8.4 Toeplitz 



(127) 

(128) 
(129) 
(130) 



(131) 



Like symmetric matrices and diagonal matrices also Toeplitz matrices has a 
special structure which should be taken into account when the derivative with 
respect to a matrix with Toeplitz structure. 



<9Tr(AT) 

<9T 
dTr(TA) 

dT 

Tr(A) Tr([A T ] rel ) Tr ( [ [A T ] 1 „ ] „ _ 1 1 2 ) 

Tr([A T ]i„)) Tr(A) 
Tr([[A T ] 1 „] 2 „_ 1 ) 



Tr([[A J ] 1?l ] 2 j) Tr([A']l„)) 



a (A) 



(132) 



Tr([[A T ] 1 „]„_ 12 ) 

Tr([A T ] rel ) 
Tr(A) 



As it can be seen, the derivative ct(A) also has a Toeplitz structure. Each value 
in the diagonal is the sum of all the diagonal valued in A, the values in the 
diagonals next to the main diagonal equal the sum of the diagonal next to the 
main diagonal in A T . This result is only valid for the unconstrained Toeplitz 
matrix. If the Toeplitz matrix also is symmetric, the same derivative yields 



dTr(AT) _ dTr(TA) 



= a(A) + a(A) T -a(A)oI 



(133) 
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3 Inverses 



3.1 Basic 

3.1.1 Definition 

The inverse A -1 of a matrix A £ 



is defined such that 



AA 1 = A _1 A = I, 



(134) 

where I is the nxn identity matrix. If A -1 exists, A is said to be nonsingular. 
Otherwise, A is said to be singular (see e.g. |HJ). 

3.1.2 Cofactors and Adjoint 

The submatrix of a matrix A, denoted by [A]y is a (n — 1) x (n — 1) matrix 
obtained by deleting the ith row and the jth column of A. The cofactor 
of a matrix is defined as 



cof(A,i,j) = (-l)^det([A] ii ) ) 
The matrix of cofactors can be created from the cofactors 

"cof(A,l,l) ••• cof(A,l,ri)' 

cof ( A )= : cof(A,i,j) : 

cof(A,n,l) ••• cof(A,n, n) 

The adjoint matrix is the transpose of the cofactor matrix 

adj(A) = (cof(A)) T , 

3.1.3 Determinant 

The determinant of a matrix A S C nx ™ j s defined as (see [T^] ) 



det(A) = ^-lJ^ydetQAly) 



j'=i 

n 

= X)4y«rf(A,l,i). 

3.1.4 Construction 

The inverse matrix can be constructed, using the adjoint matrix, by 



1 



det(A) 

For the case of 2 x 2 matrices, see section! 1.2 



adj(A) 



(135) 



(136) 



(137) 



(138) 
(139) 



(140) 
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3 INVERSES 



3.1.5 Condition number 

The condition number of a matrix c(A) is the ratio between the largest and the 
smallest singular value of a matrix (see Section 5.3 on singular values), 

c(A) = ^± (141) 

The condition number can be used to measure how singular a matrix is. If the 
condition number is large, it indicates that the matrix is nearly singular. The 
condition number can also be estimated from the matrix norms. Here 

c(A) = ||A|| • HA" 1 !!, (142) 

where || • || is a norm such as e.g the 1-norm, the 2-norm, the oo-norm or the 
Frobenius norm (see Sec 10.4 for more on matrix norms). 

The 2-norm of A equals y/ (max(eig(A ff A))) jTSJ, p. 57]. For a symmetric 
matrix, this reduces to ||A||2 = max(|eig(A)|) [T2J p. 394]. If the matrix ia 
symmetric and positive definite, ||A||2 = max(eig(A)). The condition number 
based on the 2-norm thus reduces to 

||A|| 2 ||A- 1 || 2 = max(eig(A))max(eig(A- 1 ))^ maxfc ^ (143) 

mm(eig(A)) 

3.2 Exact Relations 

3.2.1 Basic 

(AB)- 1 = B^A 1 (144) 

3.2.2 The Woodbury identity 

The Woodbury identity comes in many variants. The latter of the two can be 
found in [12] 

(A + CBC T )~ 1 = A -1 — A _1 C(B _1 + C T A _1 C) _1 C T A _1 (145) 
(A + UBV)- 1 = A -1 - A _1 U(B _1 + VA _1 U) _:l VA _1 (146) 

If P, R are positive definite, then (see [3D]) 

(P- 1 + B T R X B) 1 B T R 1 = PB T (BPB T + R)" 1 (147) 

3.2.3 The Kailath Variant 

(A + BC)- 1 = A" 1 - A _1 B(I + CA^B^CA" 1 (148) 
See H page 153]. 
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3.2.4 The Searle Set of Identities 

The following set of identities, can be found in [25J, page 151], 



(I + A- 1 )- 1 


= A(A + I)- X 


(149) 


(A + BB T ) _1 B 


= A" 1 B(I + B T A- 1 B)- 1 


(150) 


(A^ + B- 1 )- 1 


= A(A + B) _1 B = B(A + B) _1 A 


(151) 


A- A(A + B) _1 A 


= B-B(A + B) X B 


(152) 


A _1 + B 1 


= A" 1 (A + B)B- 1 


(153) 


(I + AB)- 1 


= I - A(I + BA) _1 B 


(154) 


(1 + AB) _1 A 


= A(IiBA)" 1 


(155) 



3.2.5 Rank-1 update of Moore-Penrose Inverse 

The following is a rank-1 update for the Moore-Penrose pseudo-inverse of real 
valued matrices and proof can be found in |18j . The matrix G is defined below: 



(A + cd T )+ = A+ + G (156) 



Using the the notation 



13 = l + d T A+c (157) 

v = A+c (158) 

n = (A+) T d (159) 

w = (I-AA+)c (160) 

m = (I-A+A) T d (161) 

the solution is given as six different cases, depending on the entities ||w||, 
||m||, and (3. Please note, that for any (column) vector v it holds that v + = 
v T (v T v) _1 = ||Li 2 . The solution is: 



Case 1 of 6: If ||w|| ^ and ||m|| ^ 0. Then 



G = -vw + -(m + )V + /?(m + )V (162) 
= --nAr^vw T - y— \^mn T + T — M f, M9 mw T (163) 

ll w ll ll m ll ll m lnl w ll 

Case 2 of 6: If ||w|| = and ||m|| ^ and (3 = 0. Then 

G = -vv+A+ - (m+) T n T (164) 



i+) T n T 

1 , . ... _J. 

|m 



T A+ - — mn T (165) 
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Case 3 of 6: If ||w|| = and (3^0. Then 

1 — ,tk + 3 /HKn ... ... . \ i ii'".. , a ,•/, 



G = -mv J A + - .. ..„.. — — r ^J—^-m + v '' (A+Tv + n 

||v|| 2 ||m|| 2 + |/3| 2 V P J V P 

(166) 

Case 4 of 6: If ||w|| ^ and ||m|| = and p = 0. Then 

G = -A+nn+-vw+ (167) 
= --^A+nn^T-^vw 1 (168) 



Case 5 of 6: If ||m|| = and p ^ 0. Then 
G=iA+nw T - „ „ 2 „ f, 2 |fl|2 f^A+n + v) (K W + 1 Y (i 69 ) 



Case 6 of 6: If ||w|| = and ||m|| = and /3 = 0. Then 

G = -W+A+ - A+nn+ +v+A+nvn + (170) 

= - n— rw vvTA+ - n— rw A+nnT + | V ||^"| " vn T (171) 
v n 2 "• 1 : 



3.3 Implication on Inverses 

If (A + B)- 1 = A" 1 then AB *A = BA *B (172) 

See [25]. 

3.3.1 A PosDef identity 

Assume P, R to be positive definite and invertible, then 

(P- 1 + B T R X B) 1 B T R 1 = PB T (BPB T + R)" 1 (173) 

See [5U]. 

3.4 Approximations 

The following is a Taylor expansion if A" — > when n — > oo, 

(1 + A)" 1 =1- A + A 2 - A 3 + ... (174) 
Note the following variant can be useful 

(I + A)" 1 = -\-I--A+ \a 2 - \a 3 + ...) (175) 
c [c c c z & J 

The following approximation is from [25] and holds when A large and symmetric 

A - A (I + A) _1 A = I — A -1 (176) 

If <7 2 is small compared to Q and M then 

(Q + ^M)- 1 = Q 1 -c^Q^MQ 1 (177) 
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3.5 Generalized Inverse 
3.5.1 Definition 

A generalized inverse matrix of the matrix A is any matrix A - such that (see 

BSD 

AA A = A (178) 

The matrix A~ is not unique. 

3.6 Pseudo Inverse 
3.6.1 Definition 

The pseudo inverse (or Moore-Penrose inverse) of a matrix A is the matrix A + 
that fulfils 

I AA+A = A 
II A+AA+ = A+ 

III AA + symmetric 

IV A + A symmetric 

The matrix A + is unique and does always exist. Note that in case of com- 
plex matrices, the symmetric condition is substituted by a condition of being 
Hcrmitian. 
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3.6.2 Properties 

Assume A + to be the pseudo-inverse of A, then (See [3] for some of them) 



(A+)+ = 

(A T ) + = 

(A H r = 

(A*)+ = 

(A+A)A H = 

(A+A)A T ^ 

(cA)+ = 

A+ = 

A+ = 

(A T A)+ = 

(AA T )+ = 

A+ = 

A+ = 

(A H A)+ = 

(AA H )+ = 

(AB)+ = 

/(A H A)-/(0)I = 

/(AA ff )-/(0)I - 

where A e C" xm . 

Assume A to have full rank, then 

(AA 

(A+. 



A (179) 

(A+) T (180) 

(A+) H (181) 

(A+y (182) 

A H (183) 

A T (184) 

(l/c)A+ (185) 

(A T A)+A T (186) 

A T (AA T )+ (187) 

A+(A T ) + (188) 

(A T )+A+ (189) 

(A H A)+A H (190) 

A H (AA H )+ (191) 

A+(A ff )+ (192) 

(A H )+A+ (193) 

(A+AB) + (ABB + )+ (194) 

A+[/(AA H )-/(0)I]A (195) 

A[/(A H A) - /(0)I]A+ (196) 



+ )(AA+) = AA+ (197) 

A)(A+A) = A+A (198) 

Tr(AA+) = rank(AA+) (See [26]) (199) 

Tr(A+A) = rank(A+A) (See [35]) (200) 

(AB)+ = (A + AB) + (ABB + )+ (201) 

(A®B)+ = A+®B+ (202) 



For two matrices it hold that 



3.6.3 Construction 

Assume that A has full rank, then 

A n x n Square rank(A) = n => A + = A -1 

A n x m Broad rank(A) = n => A + = A T (AA T ) _1 

A n x m Tall rank(A) = m A + = (A T A) _1 A T 
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Assume A does not have full rank, i.e. A is n x m and rank(A) = r < min(n, to) . 
The pseudo inverse A + can be constructed from the singular value decomposi- 
tion A = UDV T , by 

A+ = V.D^U^ (203) 

where U r , D r , and V r are the matrices with the degenerated rows and columns 
deleted. A different way is this: There do always exist two matrices C n x r 
and D r x to of rank r, such that A = CD. Using these matrices it holds that 

A+ = D T (DD T )" 1 (C T C)- 1 C T (204) 

See 0. 
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4 Complex Matrices 
4.1 Complex Derivatives 

In order to differentiate an expression f(z) with respect to a complex z, the 
Cauchy-Riemann equations have to be satisfied ([7]): 

df{z) _ 0%(/(z)) , . 

and 

= . gg(/(g)) + ggC/M> (206) 

or in a more compact form: 

df(z) .df(z) 

w = , w (207) 

A complex function that satisfies the Cauchy-Riemann equations for points in a 
region R is said yo be analytic in this region R. In general, expressions involving 
complex conjugate or conjugate transpose do not satisfy the Cauchy-Riemann 
equations. In order to avoid this problem, a more generalized definition of 
complex derivative is used ([24|, [6]): 

• Generalized Complex Derivative: 

df(z) Ifdfiz) df(z)\ 

^ i =2(w- i w)' (2 ° 8) 

• Conjugate Complex Derivative 

df(z) l,df(z) df(z)\ 

^ = 2llfe +1 fe} (2 ° 9) 

The Generalized Complex Derivative equals the normal derivative, when / is an 
analytic function. For a non-analytic function such as f(z) — z* , the derivative 
equals zero. The Conjugate Complex Derivative equals zero, when / is an 
analytic function. The Conjugate Complex Derivative has e.g been used by [H] 
when deriving a complex gradient. 
Notice: 

df(z) df(z) df(z) 

"dT^W + 'W (210) 

• Complex Gradient Vector: If / is a real function of a complex vector z, 
then the complex gradient vector is given by ([HI P- 798]) 

V/W - 2^ (211) 

df(z) df(z) 
dftz 93?z ' 
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4 COMPLEX MATRICES 



• Complex Gradient Matrix: If / is a real function of a complex matrix Z, 
then the complex gradient matrix is given by ([2]) 

V/(Z) = (212) 

df(Z) df(Z) 
dWl 1 d3?Z ' 

These expressions can be used for gradient descent algorithms. 
4.1.1 The Chain Rule for complex numbers 

The chain rule is a little more complicated when the function of a complex 
u = f(x) is non-analytic. For a non-analytic function, the following chain rule 
can be applied ([7]) 

dg(u) dg du dg du* 

dx du dx du* dx 

dg du /dg*\*du* 
du dx \ du J dx 

Notice, if the function is analytic, the second term reduces to zero, and the func- 
tion is reduced to the normal well-known chain rule. For the matrix derivative 
of a scalar function <?(U), the chain rule can be written the following way: 



%(U) _ Tr((gfga) r 0U) ^ Tr((%P)W) 



dX dX dX 

4.1.2 Complex Derivatives of Traces 



(214) 



If the derivatives involve complex numbers, the conjugate transpose is often in- 
volved. The most useful way to show complex derivative is to show the derivative 
with respect to the real and the imaginary part separately. An easy example is: 

0Tr(X*) = dTr(X H ) = 

dftX <99?X 1 1 

dTrCX*} . dTv(X H ) 

l ^x- = l ^x- = 1 (216) 



Since the two results have the same sign, the conjugate complex derivative (209 1 
should be used. 



<9Tr(X) _ <9Tr(X T ) 
dftX ~ dftX 



(217) 



dTrpq . dTr(X T ) 

l ^x =l ^x- = 1 (218) 

Here, the two results have different signs, and the generalized complex derivative 



(2081 should be used. Hereby, it can be seen that (92 1 holds even if X is a 
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complex number. 



dTr(AX H ) 

MX 
,<9Tr(AX ff ) 

d3X 



<9Tr(AX*) 

Mx 

,<9Tr(AX*) 
1 <9QX 



A 
A 



= A 1 



dTr(XX H ) _ <9Tr(X ff X) 
<95iX Mx 
.dTr(XX H ) .dTr(X H X) 



25RX 
i23X 



<93X <93X 
By inserting p23| and ( |224| in p08| and ( p09| , it can be seen that 

9Tr(XX ff ) 

ax 

<9Tr(XX ff ) 



ax* 



= X* 
X 



(219) 
(220) 

(221) 
(222) 

(223) 
(224) 

(225) 
(226) 



Since the function Tr(XX^) is a real function of the complex matrix X, the 



complex gradient matrix (212) is given by 



VTr(XX ff ) = 2 



<9Tr(XX 

ax 5 



= 2X 



(227) 



4.1.3 Complex Derivative Involving Determinants 



Here, a calculation example is provided. The objective is to find the derivative of 
det(X^AX) with respect to X G £mxn_ The derivative is found with respect to 
the real part and the imaginary part of X, by use of (37 1 and (33 1, det(X ff AX) 
can be calculated as (see App. B.1.2| for details) 



<9det(X g AX) 
dX 



1 /<9det(X ff AX) ,<9det(X ff AX) 



<95RX 



<93X 



= dct(X ff AX)((X ff AX)- 1 X H A) T 

and the complex conjugate derivative yields 

<9det(X ff AX) _ 1 /<9det(X ff AX) .<9dct(X ff AX) 
9X 5 ~ 2V Mx +l d^X 

= det(X ff AX)AX(X ff AX)" 1 



(228) 



(229) 
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4.2 Higher order and non-linear derivatives 

d (Ax) H (Ax) d x ff A ff Ax 



(230) 



<9x (Bx) ff (Bx) <9x x fl B H Bx 

A H Ax x ff A H AxB H Bx 
~ x ff BBx (x^B»Bx) 2 (3) 

4.3 Inverse of complex sum 

Given real matrices A,B find the inverse of the complex sum A + zB. Form 
the auxiliary matrices 

E = A + tB (232) 
F = B — tA, (233) 

and find a value of t such that E _1 exists. Then 

(A + iB)- 1 = (1 - it)(E + iF)^ 1 (234) 
= (1 - it){(E + FE^F)- 1 - i(E + FE- 1 F)- 1 FE" 1 )(235) 
= (1 -^)(E + FE" 1 F)" 1 (I-iFE" 1 ) (236) 
= (E + FE _1 F) _1 ((I — iFE -1 ) - i(tl + FE -1 )) (237) 
= (E + FE" 1 F)" 1 (I-<FE^ 1 ) 

-i(E + FE" 1 F)" 1 (a + FE" 1 ) (238) 
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5 Solutions and Decompositions 

5.1 Solutions to linear equations 
5.1.1 Simple Linear Regression 

Assume we have data (x n ,y n ) for n — l,...,N and are seeking the parameters 
a, b € K such that j/j = ax^ + b. With a least squares error function, the optimal 
values for a, b can be expressed using the notation 



x = (xi, ...,x N ) J 



( yi ,...,y N ) T l = (l,...,l) T e 



p)Vxl 



and 



i? KX = x J x R xl =x T l R n = Vl 
Ryx = y T x Ryi = y T l 



a 






Rxi 


-i 


Rx.y 


b 




Rxi 


Rn 




Ryl 



as 



5.1.2 Existence in Linear Systems 

Assume A is n x m and consider the linear system 

Ax = b 

Construct the augmented matrix B = [A b] then 

Solution 



(239) 



(240) 



Condition 

rank(A) = rank(B) = m Unique solution x 

rank(A) = rank(B) < m Many solutions x 

rank(A) < rank(B) No solutions x 



5.1.3 Standard Square 

Assume A is square and invertible, then 

Ax = b => 

5.1.4 Degenerated Square 



= A *b 



(241) 



Assume A is n x n but of rank r < n. In that case, the system Ax = b is solved 
by 

x = A+b 

where A + is the pseudo-inverse of the rank-deficient matrix, constructed as 
described in section [3.6.31 
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5.1.5 Cramer's rule 

The equation 



Ax = b, (242) 



where A is square has exactly one solution x if the zth clement in x can be 
found as 

detB 

^=detA' (243) 
where B equals A, but the ith column in A has been substituted by b. 

5.1.6 Over-determined Rectangular 

Assume A to be n x to, n > to (tall) and rank(A) = to, then 

Ax = b =s> x = (A T A)" 1 A T b = A+b (244) 

that is if there exists a solution x at all! If there is no solution the following 
can be useful: 

Ax = b =*> x min = A+b (245) 

Now x m i„ is the vector x which minimizes ||Ax — b|| 2 , i.e. the vector which is 
"least wrong". The matrix A + is the pseudo-inverse of A. See [3]. 

5.1.7 Under-determined Rectangular 

Assume A is n x to and n < to ("broad") and rank(A) = n. 

Ax = b => x mm = A T (AA T ) _1 b (246) 

The equation have many solutions x. But x m i„ is the solution which minimizes 
1 1 Ax — b|| 2 and also the solution with the smallest norm ||x|| 2 . The same holds 
for a matrix version: Assume A is n x m, X is m x n and B is n x n, then 

AX = B =► X mm = A+B (247) 

The equation have many solutions X. But X m ;„ is the solution which minimizes 
||AX — B|| 2 and also the solution with the smallest norm ||X|| 2 . See [3J. 

Similar but different: Assume A is square n x n and the matrices Bo,Bi 
are n x N, where N > n, then if Bo has maximal rank 

AB = Bi A min — BiBg (B Bq ) _1 (248) 

where A min denotes the matrix which is optimal in a least square sense. An 
interpretation is that A is the linear approximation which maps the columns 
vectors of B into the columns vectors of Bi. 

5.1.8 Linear form and zeros 

Ax = 0, Vx =s> A = (249) 
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5.1.9 Square form and zeros 

If A is symmetric, then 

x T Ax = 0, Vx =>• A = (250) 

5.1.10 The Lyapunov Equation 

AX + XB = C (251) 
vec(X) = (I <g> A + B T ® I) _1 vec(C) (252) 

Sec |10.2.1| and |10.2.2| for details on the Kronecker product and the vec op- 
erator. 

5.1.11 Encapsulating Sum 

£„A n XB„ = C (253) 
vec(X) = (E„B^A„)~ 1 vec(C) (254) 

See Sec |10.2.1| and |10.2.2| for details on the Kronecker product and the vec 
operator. 

5.2 Eigenvalues and Eigenvectors 

5.2.1 Definition 

The eigenvectors v and eigenvalues A are the ones satisfying 

Avi = XiVi (255) 

AV = VD, (D) iS = 5^, (256) 
where the columns of V are the vectors 

5.2.2 General Properties 

Assume that A E R nxm and B <E R mxn , 

eig(AB) = eig(BA) (257) 
A is n x to =4> At most min(n, m) distinct Ai (258) 
rank(A) = r ^> At most r non-zero A, (259) 
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5.2.3 Symmetric 

Assume A is symmetric, then 

VV T = I (i.e. V is orthogonal) (260) 

A,; G K (i.e. X t is real) (261) 

Tr(A p ) = Y,i>% (262) 

eig(I + cA) = 1 + cX (263) 

eig(A-cI) = Xi-c (264) 

eig(A" 1 ) = Ar 1 (265) 

For a symmetric, positive matrix A, 

eig(A T A) = cig(AA T ) = eig(A) o eig(A) (266) 

5.2.4 Characteristic polynomial 

The characteristic polynomial for the matrix A is 

= dct(A-AI) (267) 
= A" -. 9l A"- 1 + 32 A"- 2 -... + (-l)"<?„ (268) 

Note that the coefficients gj for j = l,...,n are the n invariants under rotation 
of A. Thus, gj is the sum of the determinants of all the sub-matrices of A taken 
j rows and columns at a time. That is, g\ is the trace of A, and gi is the sum 
of the determinants of the n(n — l)/2 sub-matrices that can be formed from A 
by deleting all but two rows and columns, and so on - see [17] . 

5.3 Singular Value Decomposition 

Any n x m matrix A can be written as 

A = UDV T , (269) 

where 

U = eigenvectors of AA T n x n 

D = v /diag(eig(AA T )) nxm (270) 

V = eigenvectors of A T A m x m 

5.3.1 Symmetric Square decomposed into squares 

Assume A to be n x n and symmetric. Then 

[ A ] = [ V ] [ D ] [ V r ] , (271) 

where D is diagonal with the eigenvalues of A, and V is orthogonal and the 
eigenvectors of A. 
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5.3.2 Square decomposed into squares 

Assume A e R nxn . Then 

[ A ] = [ V ] [ D ] [ & 



(272) 



where D is diagonal with the square root of the eigenvalues of AA T , V is the 
eigenvectors of AA T and U T is the eigenvectors of A T A. 

5.3.3 Square decomposed into rectangular 

Assume V*D*Uf = then we can expand the SVD of A into 



[ A ] = [ V | V, 
where the SVD of A is A = VDU T . 



D 











If 



(273) 



5.3.4 Rectangular decomposition I 

Assume A is n x to, V is n x n, D is n x n, U T is n x to 

[ A ] = [V][D][ IF ] 



(274) 



where D is diagonal with the square root of the eigenvalues of AA T , V is the 
eigenvectors of AA T and U T is the eigenvectors of A T A. 

5.3.5 Rectangular decomposition II 

Assume A is n x to, V is n x m, D is m x m, U T is to x to 



[ A 



[ v ] 



D 



(275) 



5.3.6 Rectangular decomposition III 

Assume A is n x to, V is n x n, D is n x m, U T is to x to 



[ A ] = [V][ D ] 



(276) 



where D is diagonal with the square root of the eigenvalues of AA T , V is the 
eigenvectors of AA T and U T is the eigenvectors of A T A. 
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5.4 Triangular Decomposition 

5.5 LU decomposition 

Assume A is a square matrix with non-zero leading principal minors, then 

A = LU (277) 

where L is a unique unit lower triangular matrix and U is a unique upper 
triangular matrix. 

5.5.1 Cholesky-decomposition 

Assume A is a symmetric positive definite square matrix, then 

A = U T U = LL T , (278) 

where U is a unique upper triangular matrix and L is a unique lower triangular 
matrix. 

5.6 LDM decomposition 

Assume A is a square matrix with non-zero leading principal minor^] then 

A = LDM T (279) 

where L, M are unique unit lower triangular matrices and D is a unique diagonal 
matrix. 

5.7 LDL decompositions 

The LDL decomposition are special cases of the LDM decomposition. Assume 
A is a non-singular symmetric definite square matrix, then 

A = LDL T = L T DL (280) 

where L is a unit lower triangular matrix and D is a diagonal matrix. If A is 
also positive definite, then D has strictly positive diagonal entries. 

1 If the matrix that corresponds to a principal minor is a quadratic upper-left part of the 
larger matrix (i.e., it consists of matrix elements in rows and columns from 1 to k), then the 
principal minor is called a leading principal minor. For an n times n square matrix, there are 
n leading principal minors. 1311 
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6 Statistics and Probability 

6.1 Definition of Moments 

Assume x£l" xl is a random variable 

6.1.1 Mean 

The vector of means, m, is denned by 

(m)i - ( Xi ) (281) 

6.1.2 Covariance 

The matrix of covariance M is defined by 

(M) w = (( Xi - ( Xi ))( Xj - (xj))) (282) 

or alternatively as 

M= ((x-m)(x-m) T ) (283) 

6.1.3 Third moments 

The matrix of third centralized moments - in some contexts referred to as 
coskewness - is defined using the notation 

™S = ~ " (xj)){xk - (xk))) (284) 

as 



M, = 



(285) 



where ':' denotes all elements within the given index. M 3 can alternatively be 
expressed as 

M 3 = ((x-m)(x-m) T (g)(x-m) T ) (286) 
6.1.4 Fourth moments 

The matrix of fourth centralized moments - in some contexts referred to as 
cokurtosis - is defined using the notation 



in 



(4) 
"ijkl 

as 



((xi - (a; i ))(a; j - (xj)){x k - (x k })(x t - (x t ))) (287) 



M 4 = 



(4) (4) (4) | (4) (4) (4) | | (4) (4) (4) 

m:.,,m:., ...m..L m:.;,m:.„...m:.;, ... m:.,„m:.,' ...m l ; 



(288) 



i,. n 111.. 21 •■•»(.::„! I liu.. 12'"-: : 22---" t ::n2l--- 1 '"::ln" l ::2n-" l ::im 

or alternatively as 

M 4 = ((x - m)(x - m) T ® (x - m) T <g> (x - m) T ) (289) 
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6.2 Expectation of Linear Combinations 
6.2.1 Linear Forms 

Assume X and x to be a matrix and a vector of random variables. Then (see 
See [IS]) 

£[AXB + C] = A£[X]B + C (290) 
Var[Ax] = AVar[x]A T (291) 
Cov[Ax,By] = ACov[x,y]B T (292) 

Assume x to be a stochastic vector with mean m, then (see [7]) 

E[Ax + b] = Am + b (293) 
E[Ax] = Am (294) 
£[x + b] = m + b (295) 



6.2.2 Quadratic Forms 

Assume A is symmetric, c = E[x\ and £ = Var[x]. Assume also that all 
coordinates Xi are independent, have the same central moments n±, /i2, /i3, fJ> 4 
and denote a = diag(A). Then (See [2"B"] 1 

£[x T Ax] = Tr(AS) + c T Ac (296) 
Var[x T Ax] = 2 ) u 2 Tr(A 2 ) + 4^ 2 c T A 2 c + 4/i 3 c T Aa + (fj, 4 - 3^)a T a (297) 

Also, assume x to be a stochastic vector with mean m, and covariance M. Then 
(see 0) 

£[(Ax + a)(Bx + b) T ] = AMB T + (Am + a)(Bm + b) T (298) 



E[xx T ] = M + mm T (299) 

£[xa T x] = (M + mm T )a (300) 

£[x T ax T ] = a T (M + mm T ) (301) 

£[(Ax)(Ax) T ] = A(M + mm T )A T (302) 

£[(x + a)(x + a) T ] = M+ (m + a)(m + a) T (303) 



£[(Ax + a) T (Bx + b)] = Tr(AMB T ) + (Am + a) T (Bm + b) (304) 

£[x T x] = Tr(M) + m T m (305) 

£[x T Ax] = Tr(AM) + m T Am (306) 

£[(Ax) T (Ax)] = Tr(AMA T ) + (Am) T (Am) (307) 

£[(x + a) T (x + a)] = Tr(M) + (m + a) T (m + a) (308) 



See |7]. 
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6 STATISTICS AND PROBABILITY 



6.2.3 Cubic Forms 

Assume x to be a stochastic vector with independent coordinates, mean m, 
covariance M and central moments v 3 = E[(x — m) 3 ]. Then (see [7]) 

£[(Ax + a)(Bx + b) T (Cx + c)] = Adiag(B T C)v 3 

+Tr(BMC T )(Am + a) 
+AMC T (Bm + b) 

+(AMB T + (Am + a)(Bm + b) T )(Cm + c) 
£[xx T x] = v 3 + 2Mm + (Tr(M) + m T m)m 
£[(Ax + a)(Ax + a) T (Ax + a)] = Adiag(A T A)v 3 

+ [2AMA T + (Ax + a) (Ax + a) T ](Am + a) 
+Tr(AMA T )(Am + a) 
£[(Ax + a)b T (Cx + c)(Dx + d) T ] = (Ax + a)b T (CMD T + (Cm + c)(Dm + d) T ) 

+(AMC T + (Am + a)(Cm + c) T )b(Dm + d) T 
+b T (Cm + c)(AMD T - (Am + a)(Dm + d) T ) 

6.3 Weighted Scalar Variable 

Assume x S M" xl is a random variable, w S M. nxl is a vector of constants and 
y is the linear combination y = w T x. Assume further that m, M2, M 3 , M4 
denotes the mean, covariance, and central third and fourth moment matrix of 
the variable x. Then it holds that 
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((y 
{(y 
((y 



(y) 
(y)) 2 ) 
(y)) 3 ) 
(y)) 4 ) 



T 

w m 

w T M 2 w 
w T M 3 w £g) w 

W T M4W ® w (g) w 



(309) 
(310) 
(311) 
(312) 



7 MULTIVARIATE DISTRIBUTIONS 



7 Multivariate Distributions 

7.1 Cauchy 

The density function for a Cauchy distributed vector t € R Pxl , is given by 

m,,^^ T M . . d ! t(s) -' /a „ w)/2 (313) 



r (V2) [l + ft-ziJTs-^t-M)] 



where \i is the location, S is positive definite, and T denotes the gamma func- 
tion. The Cauchy distribution is a special case of the Student-t distribution. 

7.2 Dirichlet 

The Dirichlet distribution is a kind of "inverse" distribution compared to the 
multinomial distribution on the bounded continuous variate x = [x\, . . . , xp\ 

ma p. 44] 

Up r K) P 

7.3 Normal 

The normal distribution is also known as a Gaussian distribution. See sec. [8] 

7.4 Normal-Inverse Gamma 

7.5 Gaussian 

See sec. |U 

7.6 Multinomial 

If the vector n contains counts, i.e. (n)j € 0, 1, 2, then the discrete multino- 
mial disitrbution for n is given by 

P(n|a,n) = - T ^— {[a?, J> = n (314) 
n\. . . . 

where a% are probabilities, i.e. < a,; < 1 and <Zj = 1. 

7.7 Student's t 



The density of a Student-t distributed vector t £ M Pxl , is given by 

I>/2) [i + „-i (t - / ,)TS- 1 (t - m)]^' 1 



Kt|M,s^) = M- F/2 ?SS- det(S) ? /2 (315) 
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7 MULTIVARIATE DISTRIBUTIONS 



where is the location, the scale matrix £ is symmetric, positive definite, v 
is the degrees of freedom, and T denotes the gamma function. For v = 1, the 



Student-t distribution becomes the Cauchy distribution (see sec 7.11. 

7.7.1 Mean 

E(t) = /x, v > 1 (316) 

7.7.2 Variance 

cov(t) = — v > 2 (317) 
v — 2 

7.7.3 Mode 

The notion mode meaning the position of the most probable value 

mode(t) = n (318) 

7.7.4 Full Matrix Version 

If instead of a vector t £ M Pxl one has a matrix T € R PxJV , then the Studcnt-t 
distribution for T is 

^det(n)-" /2 det(S)- Ar/2 x 

det [ft" 1 + (T - M)S" 1 (T - M) T ] ~ (,y+P)/2 (319) 

where M is the location, f2 is the rescaling matrix, S is positive definite, v is 
the degrees of freedom, and T denotes the gamma function. 

7.8 Wishart 

The central Wishart distribution for M € M PxP , M is positive definite, where 
m can be regarded as a degree of freedom parameter [TBI equation 3.8.1] [8, 
section 2.5], [IT] 



1 

2 mP/2 7r P(P-l)/4 rjjr[i(m+l -p)} 

det(E)- m / 2 det(M)( m - p - 1 )/ 2 x 



exp 



-^Tr^^M) 



(320) 



7.8.1 Mean 

E(M) = m£ (321) 
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7.9 Wishart, Inverse 

The (normal) Inverse Wishart distribution for M S ]R PxP , M is positive defi- 
nite, where m can be regarded as a degree of freedom parameter 



p(M|X3,m) 



1 



2 mP/2 7r p(p-i)/ 4 j]p r[i(m + i -p)} 

det(S) m / 2 det(M)^(" l - p - 1 )/ 2 x 
I, 



exp 



-Ti^SM" 1 ) 



(322) 



7.9.1 Mean 



E(M) = £ 



m — P — 1 



(323) 
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8 Gauss ians 



8.1 Basics 

8.1.1 Density and normalization 

The density of x <~ Af(m, S) is 



p(x) 



^det(27rS) 



exp 



--(x-mfS-'lx-m) 



(324) 



Note that if x is d-dimensional, then det(27rS) = (27r) d det(S). 
Integration and normalization 



/ 



/ 



exp 
exp 



(x-m) T S- 1 (x-m) 



/ 



x'E^x + m'E^x 
1 



2 

exp 



dx = ^det(27rS) 
dx = ^/det(2-7rS) exp 



— -x T Ax + c T x 
2 



dx = ^/det(2-7rA _1 ) exp 



If X = [xix 2 ...x„] and C = [cic 2 ...c„], then 



/ 



exp 



--Tr(X T AX) + Tr(C T X) 



rfX = ^det(27rA- 1 )" exp 



c T A~ T c 



-Tr^A^C 



The derivatives of the density are 
dp(x) 



-/>(x)E _1 (x - m) 
. )x . )xJ P(x) (S-^X - m)(x - m)^- 1 - E" 1 ) 



<9x 

d 2 p 



(325) 
(326) 



8.1.2 Marginal Distribution 

Assume x ~ 7V x (/i, E) where 



x b 



A*6 



S a E c 



then 



p(x a ) = 7V Xa (/x a ,S a ) 
p(x b ) = 7V X6 (/x h ,S 6 ) 



(327) 



(328) 
(329) 
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8 GAUSSIANS 



8.1.3 Conditional Distribution 

Assume x ~ A/" x (/i, S) where 



^6 



then 



p(x a |x h ) = W Xa (/ta,£ a ) 
p(x fc |x a ) = Af Xb (p. b , E b ) 



r Aa 
1 s» 



— ^c^f, E^T 

St = Ef, — S^Sq 1 S C 



(330) 

(331) 
(332) 



Note, that the covariance matrices are the Schur complement of the block ma- 
trix, see |9.1.5| for details. 



8.1.4 Linear combination 

Assume x ~ Af(m x , E x ) and y ~ Af(m y , E a ) then 

Ax + By + c - jVfAnij + Bm y + c, AT, X A T + BE. y B T ) 



(333) 



8.1.5 Rearranging Means 

V / det(2 7 r(A T S- 1 A)- r 



AA Ax [m,S] = 



v /det(27rS) 



AA x [A- 1 m,(A T S- 1 A)- 1 ] (334) 



8.1.6 Rearranging into squared form 

If A is symmetric, then 
1 



1 



VAx+b J x 



Tr(X T AX) + Tr(B T X) 



-(x - A-Md) t A(x - A- J b) + -b T A *b 

2 V ' y ' 2 

^Tr[(X - A- X B) T A(X - A^B)] + ^Tr(B T A" 



8.1.7 Sum of two squared forms 

In vector formulation (assuming Ei,E2 are symmetric) 

- -V 



m 1 ) T Sr 1 (x- 


mi) 


(335) 


m 2 ) T S 2 - 1 (x- 


m 2 ) 


(336) 




m c ) + C* 


(337) 
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S- 1 = Er'+Ea 1 
m c = (S^ 1 + S^ 1 )- 1 (S ] ; 1 m 1 + S^ 1 m 2 ) 



(338) 
(339) 



C = -(mfSr 1 +m^S^ 1 )(Sr 1 + S 2 - 1 )- 1 (Sr 1 m 1 + S 2 - 1 m 2 )(340) 



i(mfS 1 1 m 1 + m 2 fS 2 1 m 2 ) 



In a trace formulation (assuming Si, S 2 are symmetric) 

1. 



Tr((X - Mi) J Sr (X - Mi)) 



Tr((X-M 2 ) T S 2 - 1 (X-M 2 )) 



= --Tr[(X-M c ) T S- 1 (X-M c )] + C 



(341) 

(342) 
(343) 
(344) 



sr 1 



sr 1 



1 ^ ^2 

m c = (sr 1 + s 2 - 1 )- 1 (sr 1 Mi + s 2 ; 1 M 2 ) 
1 



(345) 
(346) 



C 



-TV 



1 



(Sr'Mi + S 2 r 1 M 2 ) T (S^ 1 + SJ 1 )- 1 ^^! + E^Ma) 

' V ! A I T>, I ' \ ■ A r i (347) 



--Tr(Mf £r Mi + ^ E 2 'M,) 



8.1.8 Product of gaussian densities 

Let A/" x (m, S) denote a density of x, then 

A/x(mi,Ei) ■A^K.Sj) = c c 7V x (m c ,S c ) 



(348) 



JV mi (m2,(£i + £ 2 )) 
1 



: exp 



--(mi - m 2 ) T (Si + S 2 ) 1 (nii-m 2 ) 



v /det(27r(S 1 + S 2 )) 
m c = (E^ 1 + S 2 " 1 )- 1 (S 1 ; 1 mi + E 2 T 1 m 2 ) 

s c = (e^ + e,- 1 )- 1 

but note that the product is not normalized as a density of x. 
8.2 Moments 

8.2.1 Mean and covariance of linear forms 

First and second moments. Assume x <~ 7V(m, E) 

£(x) = m 



(349) 
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Cov(x, x) = Var(x) = S = £(xx T ) - E(x)E(x T ) = E(xx T ) - mm T (350) 

As for any other distribution is holds for gaussians that 

B[Ax] = AE[x] (351) 
Var[Ax] = AVar[x]A T (352) 
Cov[Ax,By] = ACov[x,y]B T (353) 

8.2.2 Mean and variance of square forms 

Mean and variance of square forms: Assume x ~ A/"(m, E) 

E(xx T ) = S + mm T (354) 

£[x T Ax] = Tr(AS) + m T Am (355) 
Var(x T Ax) = Tr[AE(A + A T )E] + ... 

+m T (A + A T )£(A + A T )m (356) 

£[(x - m') T A(x - m')] = (m — m') T A(m — m') + Tr(AE) (357) 

If E = cr 2 I and A is symmetric, then 

Var(x T Ax) = 2cr 4 Tr(A 2 ) + 4cr 2 m T A 2 m (358) 

Assume x <~ Af(0, cr 2 I) and A and B to be symmetric, then 

Cov(x T Ax,x T Bx) = 2cr 4 Tr(AB) (359) 

8.2.3 Cubic forms 

E [xb T xx T ] = mb T (M + mm 7 ' ) + (M + mm T ) bm T 

+b T m(M - mm T ) (360) 

8.2.4 Mean of Quartic Forms 

E[xx T xx T ] = 2(S + mm T ) 2 + m T m(S - mm T ) 
+Tr(S)(E + mm T ) 
£[xx T Axx T ] = (S + mm T )(A + A T )(£ + mm T ) 

+m T Am(S - mm T ) + Tr[AS](S + mm T ) 
E[x T xx T x] = 2Tr(S 2 ) + 4m T Sm + (Tr(S) + m T m) 2 
£[x T Axx r Bx] = Tr[AS(B + B T )£] +m T (A + A T )S(B + B T )m 
+(Tr(AS) + m T Am)(Tr(BS) + m T Bm) 
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J B[a T xb T xc T xd T x] 
= (a T (52 + mm T )b)(c T (52 + mm T )d) 
+ (a T (52 + mm T )c)(b T (52 + mm T )d) 

+ (a T (52 + mm T )d)(b T (E + mm T )c) - 2a T mb 7 mc 7 md 7 m 

E[(Ax + a)(Bx + b) T (Cx + c)(Dx + d) T ] 
= [A52B T + (Am + a)(Bm + b) T ][C52D T + (Cm + c)(Dm + d) T ] 
+ [A52C T + (Am + a)(Cm + c) T ][B52D T + (Bm + b)(Dm + d) T ] 
+ (Bm + b) T (Cm + c) [ASD T - (Am + a)(Dm + d) T ] 
+Tr(B52C T )[A52D T + (Am + a) (Dm + d) T ] 

£[(Ax + a) T (Bx + b)(Cx + c) T (Dx + d)] 
= Tr[A52(C T D + D T C)£B T ] 

+ [(Am + a) T B + (Bm + b) T A]52[C T (Dm + d) + D T (Cm + c)] 

+ [Tr(A52B T ) + (Am + a) T (Bm + b)] [Tr(C52D T ) + (Cm + c) T (Dm + d)] 

See |7]. 

8.2.5 Moments 

E[x] = 5> fc m fc (361) 

k 

Cov(x) = ^2 ^2 PkPk' (52 fe + m k ml - m^m^, ) (362) 

k k' 

8.3 Miscellaneous 

8.3.1 Whitening 

Assume x ~ A/"(m, 52) then 

z = S- x / 2 (x - m) - Af(0, 1) (363) 

Conversely having z ~ Af(0, 1) one can generate data x ~ Af(m, 52) by setting 

x = 52 1/2 z + m ~ A/"(m, 52) (364) 

Note that 52 1 / 2 means the matrix which fulfils E 172 !] 1 / 2 = 52, and that it exists 
and is unique since 52 is positive definite. 

8.3.2 The Chi-Square connection 

Assume x ~ Af(m, 52) and x to be n dimensional, then 



z 



(x-m) T E- 1 (x-m)^ X 2 (365) 



where \n denotes the Chi square distribution with n degrees of freedom. 
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8.3.3 Entropy 

Entropy of a D-dimensional gaussian 



H(x) = - J W(m, S) ln7V(m, S)rfx = In v/det(27r£) 



2~ 



(366) 



8.4 Mixture of Gaussians 
8.4.1 Density 

The variable x is distributed as a mixture of gaussians if it has the density 

K 



p(x) = p k — — exp 



1 



(x-m fc ) T S fc 1 (x-m fe ) 



(367) 



where pfc sum to 1 and the all are positive definite. 
8.4.2 Derivatives 

Defining p(s) = J2k PkN s {Vk> s fc) one get 
ainp(s) pjK(fij,^j) d 



d Pj 



dlnp(s) 



d\np(s) 
<95L 



p j Af s (n j ,'E J ) 1 

EfePfeA/' s (^fe,Sfe) Pj 
PjK(j^£j) d 



pjMjjj^j) d 



ln^-TV; (/i,, E,)] 



(368) 
(369) 
(370) 
(371) 
(372) 



E fe p fc A4(/x fc ,S fe )2 L ' 
But pfc and Sfc needs to be constrained 
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9 Special Matrices 

9.1 Block matrices 

Let Aij denote the ijth block of A. 

9.1.1 Multiplication 

Assuming the dimensions of the blocks matches we have 



A n 


A12 ' 


A 2 i 


A 22 . 



Bn 


B12 




' AnBu- 


I- A12B21 


A11B12 - 


- A12B22 


B 2 i 


B22 




A 2 iB n - 


- A 22 B 2 i 


A2iB 12 - 


~ A22B22 



9.1.2 The Determinant 

The determinant can be expressed as by the use of 

Ci = An — A 12 A 22 1 A 2 i 
C 2 = A 2 2 — A 2 iA 11 Ai2 



(374) 
(375) 



as 



det 



A n 


A12 ' 


A 2 i 


A 22 . 



det(A 22 ) • det(Ci) = det(An) • det(C 2 ) 



9.1.3 The Inverse 

The inverse can be expressed as by the use of 

Ci 



as 



An 


A12 ' 


A 2 i 


A 22 . 



An 


- AiaA^ 1 


A 2 i 




A 22 


- A21 A^ 1 


A12 






cr 1 


A 


11 Ai 2 C 2 1 


[ -C2 


'AaiAJi 1 


c- J 



(376) 
(377) 



An 1 - 


hA n 1 A 12 C 2 - 1 A 2 iA n 1 


—C 1 1 Ai 2 A 22 1 


-A 22 A2iCj 




hA 22 1 A 2 iC 1 1 A 12 A 22 1 _ 



9.1.4 Block diagonal 

For block diagonal matrices we have 



det 



An 





-1 


r (An)- 1 








A 22 _ 







(A22)- 1 



An 








A 22 . 



det(An) •det(A 22 ) 



(378) 
(379) 
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9.1.5 Schur complement 

The Schur complement of the matrix 



An 


A12 ' 


A 2 i 


A 22 



is the matrix 



An - Ai 2 A 22 1 A 2 i 



that is, what is denoted C 2 above. Using the Schur complement, one can rewrite 
the inverse of a block matrix 



An 


A12 ' 


A 2i 


A 22 _ 



n -1 



A,, A 2 j 



(An - Ai 2 A 22 1 A 21 ) 







-AiaA^ 1 



The Schur complement is useful when solving linear systems of the form 



An 


A12 




Xl 




r bi l 


A 2i 


A 22 _ 




. X 2 




b 2 



which has the following equation for xi 

(An - Ai 2 A 22 A 2 i)xi = bi - Ai 2 A 22 b 2 

When the appropriate inverses exists, this can be solved for Xi which can then 
be inserted in the equation for x 2 to solve for x 2 . 

9.2 Discrete Fourier Transform Matrix, The 

The DFT matrix is an TV x N symmetric matrix Wat, where the k, nth element 
is given by 

W kn = e ^^2 ( 380 ) 

Thus the discrete Fourier transform (DFT) can be expressed as 

N-l 



X(k) = J2 <n)W k N n . 



(381) 



n=0 



Likewise the inverse discrete Fourier transform (IDFT) can be expressed as 



JV-l 



(") = Jf E X{k)W^. 



(382) 



fc=0 



The DFT of the vector x = [x(0), x(l), • • • , x(N — 1)] T can be written in matrix 
form as 

X = W N x, (383) 
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where X = [X(0),X(1), • • • , x(N - 1)] T . The IDFT is similarly given as 

x = W^X. (384) 

Some properties of exist: 

W^ 1 - ^W^ (385) 

WatW^ = NI (386) 

= W$ (387) 



— j'27r 

If Wn — eT ; then 



W ™+N/2 = _ w rn (388) 



N 

Notice, the DFT matrix is a Vandermonde Matrix. 

The following important relation between the circulant matrix and the dis- 
crete Fourier transform (DFT) exists 

Tc=W w 1 (Io(W A rt))W JV , (389) 

where t = [to, t\, ■ ■ ■ , t n ^i] T is the first row of Tq- 

9.3 Hermitian Matrices and skew-Hermitian 

A matrix A e c mx ™ j s called Hermitian if 

A H = A 

For real valued matrices, Hermitian and symmetric matrices are equivalent. 

A is Hermitian <S> x"Ax G R, Vx E C" x1 (390) 
A is Hermitian <^4> eig(A) e K (391) 

Note that 

A = B + iC 

where B,C are hermitian, then 

_, A + A H „ A-A ff 
B 



2 ' 2i 
9.3.1 Skew-Hermitian 

A matrix A is called shew -hermitian if 

A = —A H 

For real valued matrices, skew-Hcrmitian and skew-symmetric matrices are 
equivalent. 

A Hermitian <^> iA is skew-hermitian (392) 
A skew-Hermitian <f> x H Ay = -x ff A H y, Vx,y (393) 
A skew-Hermitian eig(A) = iA, AeR (394) 
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9.4 Idempotent Matrices 

A matrix A is idempotent if 

AA = A 

Idempotent matrices A and B, have the following properties 

A" = A, forn= 1,2,3,... (395) 

I A is idempotent (396) 

A H is idempotent (397) 

I A H is idempotent (398) 

If AB = BA => AB is idempotent (399) 

rank(A) = Tr(A) (400) 

A(I-A) = (401) 

(I-A)A = (402) 

A+ = A (403) 

/(al + tA) = (I-A)/(s) + A/(« + i) (404) 

Note that A — I is not necessarily idempotent. 

9.4.1 Nilpotent 

A matrix A is nilpotent if 

A 2 = 

A nilpotent matrix has the following property: 

/(sI + iA) = If (a) + tAf'(s) (405) 

9.4.2 Unipotent 

A matrix A is unipotent if 

AA = I 

A unipotent matrix has the following property: 

/(sI + iA) = [(I + A)/(* + t) + (I-A)/(s-t)]/2 (406) 

9.5 Orthogonal matrices 

If a square matrix Q is orthogonal, if and only if, 

Q T Q = QQ T = I 

and then Q has the following properties 

• Its eigenvalues are placed on the unit circle. 
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• Its eigenvectors are unitary, i.e. have length one. 

• The inverse of an orthogonal matrix is orthogonal too. 
Basic properties for the orthogonal matrix Q 

Q 1 = Q T 

Q T = Q 

QQ T = I 

Q T Q = I 

det(Q) = ±1 

9.5.1 Ortho-Sym 

A matrix Q + which simultaneously is orthogonal and symmetric is called an 
ortho-sym matrix |20j . Hereby 

QTQ+ = I (407) 

Q+ = Ql (408) 

The powers of an ortho-sym matrix are given by the following rule 

= ^^I + ^±^Q + (409) 
= l + coB(fe7r) I+ l-cos(fc7r) Q+ (41Q) 

9.5.2 Ortho-Skew 

A matrix which simultaneously is orthogonal and antisymmetric is called an 
ortho-skew matrix [50]. Hereby 

Q H Q = I (411) 
Q- = Q H (412) 

The powers of an ortho-skew matrix are given by the following rule 

= cos(fc|)I + sin(fc|)Q_ (414) 

9.5.3 Decomposition 

A square matrix A can always be written as a sum of a symmetric A + and an 
antisymmetric matrix A_ 

A = A+ + A_ (415) 
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9.6 Positive Definite and Semi-definite Matrices 

9.6.1 Definitions 

A matrix A is positive definite if and only if 

x T Ax > 0, Vx ^ (416) 
A matrix A is positive semi-definite if and only if 

x T Ax > 0, Vx (417) 
Note that if A is positive definite, then A is also positive semi-definite. 

9.6.2 Eigenvalues 

The following holds with respect to the eigenvalues: 



A pos. def. ^ eig( A+ ^" ) > 
A pos. semi-def. eig( A+ 2 A ) > 

9.6.3 Trace 

The following holds with respect to the trace: 

A pos. def. Tr(A) > 

A pos. semi-def. =^> Tr(A) > 



(418) 



(419) 



9.6.4 Inverse 

If A is positive definite, then A is invertible and A -1 is also positive definite. 

9.6.5 Diagonal 

If A is positive definite, then An > 0,Vi 

9.6.6 Decomposition I 

The matrix A is positive semi-definite of rank r <=> there exists a matrix B of 
rank r such that A — BB T 

The matrix A is positive definite there exists an invertible matrix B such 
that A = BB T 

9.6.7 Decomposition II 

Assume A is an n x n positive semi-definite, then there exists an n x r matrix 
B of rank r such that B T AB = I. 
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9.6.8 Equation with zeros 

Assume A is positive semi-definite, then X T AX = =>■ AX = 

9.6.9 Rank of product 

Assume A is positive definite, then rank(BAB T ) = rank(B) 

9.6.10 Positive definite property 

If A is n x n positive definite and B is r x n of rank r, then BAB T is positive 
definite. 

9.6.11 Outer Product 

If X is n x r, where n < r and rank(X) = n, then XX T is positive definite. 

9.6.12 Small pertubations 

If A is positive definite and B is symmetric, then A — tB is positive definite for 
sufficiently small t. 

9.6.13 Hadamard inequality 

If A is a positive definite or semi-definite matrix, then 

det(A) < Y[A U 

i 

See [13 pp.477] 

9.6.14 Hadamard product relation 

Assume that P = AA T and Q = BB T are semi positive definite matrices, it 
then holds that 

P o Q = RR T 

where the columns of R are constructed as follows: r i+ ^_i^ A — a; o b,, for 
i = 1,2,..., Na and j = 1, 2, Nb- The result is unpublished, but reported by 
Pavel Sakov and Craig Bishop. 

9.7 Singleentry Matrix, The 
9.7.1 Definition 

The single-entry matrix J lJ 6 E" xn is defined as the matrix which is zero 
everywhere except in the entry in which it is 1. In a 4 x 4 example one 
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might have 



r23 





10 







(420) 



The single-entry matrix is very useful when working with derivatives of expres- 
sions involving matrices. 

9.7.2 Swap and Zeros 

Assume A to be n x m and JP- 7 to be m x p 

A3 ij = [ ... A, ... ] (421) 

i.e. an n x p matrix of zeros with the i.th column of A in place of the j.th 
column. Assume A to be n x m and J y to be p x n 





J^A 





A, 







(422) 



i.e. an p x m matrix of zeros with the j.th row of A in the placed of the i.th 
row. 

9.7.3 Rewriting product of elements 



A ki Bji = (Ae.ejB) 
A ik B tj = (A T e,ejB T ) 



u — 



AikBji = (A T e 4 ejB) fe; 



(AJ«B) fc , 

(A T 3^B T ) kl 

(A T J«B) H 



(AJ^B T ) 



9.7.4 

If i = j 



A ki B l3 = (Ae;ej B J ) kl 
Properties of the Singleentry Matrix 

.J J = .) (J lj ) T (3^) T = J 

jij (jij)T = jij (jij jTjij = jij 

jijjij = o (J^') T (J«) T = 

jij (jij')T = jii (J«) T jiJ = Jij 



kl 



(423) 
(424) 
(425) 
(426) 
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9.7.5 The Singleentry Matrix in Scalar Expressions 

Assume A is n x to and J is to x n, then 

Tr(AJ 4J ) = Tr(J^A) = (A T )„- (427) 

Assume A is n x n, J is n x to and B is to x n, then 

Tr(AJ ij 'B) = (A T B T )y (428) 
Tr(AJ J '*B) = (BA)jj (429) 
Tr(AJ i - 7 'j ij 'B) = diag(A T B T ) u - (430) 



Assume A is n x n, J IJ is n x to B is to x n, then 



x T AJ y Bx = (A T xx T B T ) y (431) 
x T AJ ij J ij 'Bx = diag(A T xx T B T ) u - (432) 



9.7.6 Structure Matrices 

The structure matrix is defined by 

H = s " < 433 > 

If A has no special structure then 

S ij = J i: > (434) 

If A is symmetric then 

S" = .) + .1 - J .1" (435) 

9.8 Symmetric, Skew-symmetric/ Antisymmetric 

9.8.1 Symmetric 

The matrix A is said to be symmetric if 

A = A T (436) 

Symmetric matrices have many important properties, e.g. that their eigenvalues 
are real and eigenvectors orthogonal. 

9.8.2 Skew-symmetric / Antisymmetric 

The antisymmetric matrix is also known as the skew symmetric matrix. It has 
the following property from which it is defined 



(437) 
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Hereby, it can be seen that the antisymmetric matrices always have a zero 
diagonal. The n x n antisymmetric matrices also have the following properties. 



det(A T ) 
-det(A) 



det(-A) = (-l) n det(A) 
det(— A) =0, if n is odd 



(438) 
(439) 



The eigenvalues of an antisymmetric matrix are placed on the imaginary axis 
and the eigenvectors are unitary. 



9.8.3 Decomposition 

A square matrix A can always be written as a sum of a symmetric A + and an 
antisymmetric matrix A 

A = A+ + A (440) 
Such a decomposition could e.g. be 



A = 



A + A 1 



A- A 1 



= A+ + A_ 



(441) 



9.9 Toeplitz Matrices 

A Toeplitz matrix T is a matrix where the elements of each diagonal is the 
same. In the n x n square case, it has the following structure: 



■■ *12 



to 

t-1 

t-(n-l) 



t- 



tn.— l 



h 
t 



(442) 



A Toeplitz matrix is persymmetric. If a matrix is persymmetric (or orthosym- 
metric), it means that the matrix is symmetric about its northeast-southwest 
diagonal (anti-diagonal) 12J. Persymmetric matrices is a larger class of matri- 
ces, since a persymmetric matrix not necessarily has a Toeplitz structure. There 
are some special cases of Toeplitz matrices. The symmetric Toeplitz matrix is 
given by: 



The circular Toeplitz matrix: 



to 
ti 

-(n-l) 
to t 



ti 
t 



(443) 



(444) 



to 
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The upper triangular Toeplitz matrix: 

to t\ ■ ■ ■ t n -\ 

'•• '•• : 

: '■■ '■■ t x 

■•• t 

and the lower triangular Toeplitz matrix: 

t a ••• 

/ i '■■ '■■ : 

: '•• '■• 

t-(n-l) ' ' ' t-i to 



(445) 



(446) 



9.9.1 Properties of Toeplitz Matrices 

The Toeplitz matrix has some computational advantages. The addition of two 
Toeplitz matrices can be done with 0(n) flops, multiplication of two Toeplitz 
matrices can be done in 0(n In n) flops. Toeplitz equation systems can be solved 
in 0(n 2 ) flops. The inverse of a positive definite Toeplitz matrix can be found 
in 0(n 2 ) flops too. The inverse of a Toeplitz matrix is persymmetric. The 
product of two lower triangular Toeplitz matrices is a Toeplitz matrix. More 
information on Toeplitz matrices and circulant matrices can be found in |13[l7]. 



9.10 Transition matrices 

A square matrix P is a transition matrix, also known as stochastic matrix or 
probability matrix, if 

o < (P)« < i, J2( p h = i 



The transition matrix usually describes the probability of moving from state i 
to j in one step and is closely related to markov processes. Transition matrices 
have the following properties 

Prob[i j in 1 step] = (P) y (447) 

Prob[i j in 2 steps] = (P% (448) 

Prob[i -> j in k steps] = (P% (449) 

If all rows are identical =>• P n = P (450) 

aP = a, a: is called invariant (451) 

where a is a so-called stationary probability vector, i.e., < ci^ < 1 and a i = 
1. 
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9.11 Units, Permutation and Shift 

9.11.1 Unit vector 

Let ej € R nxl be the ith unit vector, i.e. the vector which is zero in all entries 
except the ith at which it is 1. 

9.11.2 Rows and Columns 



i.th row of A 
j.th column of A 



e/A 



(452) 
(453) 



9.11.3 Permutations 

Let P be some permutation matrix, e.g. 



[ e 2 e x e 3 



For permutation matrices it holds that 

PP T = I 

and that 

AP = [ Ae 2 Ae! Ae 3 ] 



PA 



e T 
T 
e i 



(454) 



(455) 



(456) 



That is, the first is a matrix which has columns of A but in permuted sequence 
and the second is a matrix which has the rows of A but in the permuted se- 
quence. 



9.11.4 Translation, Shift or Lag Operators 

Let L denote the lag (or 'translation' or 'shift') operator defined on a 4 x 4 
example by 

" " 
10 
10 
10 



(457) 



i.e. a matrix of zeros with one on the sub-diagonal, (L)ij = <5i.j+i- With some 
signal x t for t = 1, ...,N, the n.th power of the lag operator shifts the indices, 
i.e. 

for t = 1, .., n 
x t - n for t = n+ 1, N 



(L"x) t = { 



(458) 
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A related but slightly different matrix is the 'recurrent shifted' operator defined 
on a 4x4 example by 

" 1 " 
10 
10 
10 



(459) 



i.e. a matrix defined by (L)jj = <5ij+i + <5i,i<5j.dim(L)- On a signal x it has the 
effect 

{t n x) t = x t >, t' = [(t-ri) mod TV] + 1 (460) 

That is, L is like the shift operator L except that it 'wraps' the signal as if it 
was periodic and shifted (substituting the zeros with the rear end of the signal) . 
Note that L is invertible and orthogonal, i.e. 

I" 1 = L T (461) 



9.12 Vandermonde Matrices 

A Vandermonde matrix has the form [15) 



V = 



1 vi v{ 
1 v 2 v\ 



1 v n vl 



..n-l 



,,n-l 



(462) 



The transpose of V is also said to a Vandermonde matrix. The determinant is 
given by [29] 

detV = '[[(v l -v j ) (463) 
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10.1 Functions and Series 

10.1.1 Finite Series 

(X™ - I)(X - I)" 1 = I + X + X 2 + ... + X"- 1 (464) 

10.1.2 Taylor Expansion of Scalar Function 

Consider some scalar function /(x) which takes the vector x as an argument. 
This we can Taylor expand around xo 

/(x) S /(x ) + g(x ) T (x - xo) + * (x - x ) T H(x )(x - x ) (465) 

where 

9 /W N _ 92 /(x) 



(xo) = ^ H(x ) = 

ox 



dxdx 



T 



10.1.3 Matrix Functions by Infinite Series 

As for analytical functions in one dimension, one can define a matrix function 
for square matrices X by an infinite series 

oo 

f (X) = C « X " ( 466 ) 

71=0 

assuming the limit exists and is finite. If the coefficients c„ fulfils '}2 in c n x n < oo, 
then one can prove that the above series exists and is finite, see pQ. Thus for 
any analytical function f(x) there exists a corresponding matrix function f (x) 
constructed by the Taylor expansion. Using this one can prove the following 
results: 

1) A matrix A is a zero of its own characteristic polynomium pQ: 

p(A) = det(IA - A) = C » A " K A ) = ( 467 ) 

n 

2) If A is square it holds that [T] 

A = UBLT 1 =>■ f(A) = Uf(B)U" 1 (468) 

3) A useful fact when using power series is that 

A" -> Oforn -> oo if |A| < 1 (469) 



Petersen & Pedersen, The Matrix Cookbook, Version: November 14, 2008, Page 58 



10.2 Kronecker and Vec Operator 10 FUNCTIONS AND OPERATORS 



10.1.4 Exponential Matrix Function 



In analogy to the ordinary scalar exponential function, one can define exponen- 
tial and logarithmic matrix functions: 



-A 



oo 

V-A n = I + A+-A 2 + ... 
^ n\ 2 

oo 

S^(-1) W A» = I-A + 5 A 

n=0 

oo 

^ - (t A) n = I + tA + -t 2 A 2 



n=0 
oo 



ln(I + A) ee 



(-iy 



-A n = A - -A 2 + -A 3 
2 3 



Some of the properties of the exponential function are pQ 

A+B if AB = BA 



e A e B 



I A 



(e A ) 
d 

—t 
dt 



e 

e- A 

Ae tA = e* A A, 



t e 



dt 



Tr(e tA ) = Tr(Ae tA ) 
det(e A ) = e Tr < A ) 



10.1.5 Trigonometric Functions 



sin(A) = y \' „ =A-^A 3 +4A 5 
v ; ^ 2n+l! 3! 5! 

n— 



n=0 



(470) 
(471) 
(472) 
(473) 



(474) 
(475) 

(476) 

(477) 
(478) 



(479) 
(480) 



10.2 Kronecker and Vec Operator 
10.2.1 The Kronecker Product 

The Kronecker product of an m x n matrix A and an r x q matrix B, is an 
mr x nq matrix, A <g> B defined as 



A (g) B 



A n B A 12 B 
A21H A22B 



A ln B 
A 2n B 



Ami's A„ l2 B ... A mn B 



(481) 
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The Kronecker product has the following properties (see [T9"] ) 



(A 



A®(B + C) 


— 


A H A C 


(482) 


A <g> B 




B (g A in general 


(483) 


A® (B <8>C) 


= 


(A<gB)<gC 


(484) 


«aA <g obB) 




a A a B (A <g B) 


(485) 


(A<gB) J 




A J (g B J 


(486) 


<g>B)(C<g>D) 




AC (g BD 


(487) 


(A OB)" 1 




A" 1 (gB" 1 


(488) 


(A<g>B)+ 




A+ (g>B+ 


(489) 


rank(A <g> B) 




rank(A)rank(B) 


(490) 


Tr(A ® B) 




Tr(A)Tr(B) = Tr(A A <8 A B ) 


(491) 


det(A® B) 




det(A) rank ( B ) det(B) rank < A ) 


(492) 


{eig(A<gB)} 




{eig(B (g> A)} if A, B are square 


(493) 


{eig(A<gB)} 




{eig(A)eig(B) T } 


(494) 






if A, B are symmetric and square 




cig(A <g B) 




eig(A) (g eig(B) 


(495) 



Where {A^} denotes the set of values Xi, that is, the values in no particular 
order or structure, and Aa denotes the diagonal matrix with the eigenvalues of 
A. 



10.2.2 The Vec Operator 



The vec-operator applied on a matrix A stacks the columns into a vector, i.e. 
for a 2 x 2 matrix 



A = 



An 
A 2 i 



A12 
A22 



vec(A) 



An 
A 2 i 
A12 
A22 



Properties of the vec-operator include (see [T9"] ) 



vec(AXB) = (B T (g A)vec(X) 

Tr(A T B) = vec(A) T vec(B) 

vec(A + B) = vcc(A) + vec(B) 

vec(aA) = a • vec(A) 



(496) 
(497) 
(498) 
(499) 
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10.3 Vector Norms 
10.3.1 Examples 



Mi 



Ei 



||x||oo = max|arj| 

i 

Further reading in e.g. [12j p. 52] 

10.4 Matrix Norms 

10.4.1 Definitions 

A matrix norm is a mapping which fulfils 

l|A|| > 

||A|| = O^A = 

||cA|| = |c|||A||, ceR 

||A + B|| < ||A|| + ||B|| 

10.4.2 Induced Norm or Operator Norm 

An induced norm is a matrix norm induced by a vector norm by 

||A||=sup{||Ax|| | ||x|| = l} 

where 1 1 • 1 1 ont the left side is the induced matrix norm, while 1 1 • 
side denotes the vector norm. For induced norms it holds that 



(500) 
(501) 
(502) 
(503) 



(504) 
(505) 
(506) 
(507) 



the following 
(508) 
| on the right 



l|I| 
||Ax| 

IABI 



= 1 

< HA| 

< IIAI 



Ml 
IBI 



for all A, x 
for all A,B 



10.4.3 Examples 

l|A||i 

l|A|| 2 
HAIL 



max 

j 



maxeig(A H A) 
ix 1 1 A: 

E I A 



( max HAxIL) 1 ^ 

l|x|| P = l 



(509) 
(510) 
(511) 



(512) 

(513) 
(514) 

(515) 
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max \Aij\ 
lfc(A)||i 



Tr(AA^) (Frobenius) 



(516) 

(517) 
(518) 



|A|| KF = ||sing(A)||i (Ky Fan) 
where sing(A) is the vector of singular values of the matrix A. 

10.4.4 Inequalities 

E. H. Rasmussen has in yet unpublished material derived and collected the 
following inequalities. They are collected in a table as below, assuming A is an 
m x n, and d = rank(A) 



|A||i 



KF 



I max 

ji m 
loo n 
| 2 y/rrm 
|f y/mn 

|kf V mnd 



1 

m 



n 



rid 



m 
m 
md 



1 

\frh 

Vd 

d 



1 



m 
n 



1 

Vd 



which are to be read as, e.g. 



lAllo < vm • IIA 



(519) 



10.4.5 Condition Number 



The 2-norm of A equals \J (max(eig(A T A))) p. 5 7]. For a symmetric, pos- 
itive definite matrix, this reduces to max(eig(A)) The condition number based 
on the 2-norm thus reduces to 



lAIMIA- 1 ^ = max(eig(A))max(eig(A^ 1 )) = 



max(eig(A)) 



min(cig(A)) 

10.5 Rank 

10.5.1 Sylvester's Inequality 

If A is m x n and B is n x r, then 

rank(A) + rank(B) — n < rank(AB) < min{rank(A), rank(B)} 

10.6 Integral Involving Dirac Delta Functions 

Assuming A to be square, then 

1 



(520) 



(521) 



p(s)<5(x — As)ds 



det(A) 



p(A-! x ) 



(522) 
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10 FUNCTIONS AND OPERATORS 



Assuming A to be "underdetermined" , i.e. "tall", then 

J^=p(A+x) if x = AA+x 

-f A T A \ " \ / 



J p(b)5(x - As)ds = | Vdet(A^A) 3 



elsewhere 



(523) 



See 0. 



10.7 Miscellaneous 

For any A it holds that 

rank(A) = rank(A T ) = rank(AA T ) = rank(A T A) (524) 

It holds that 

A is positive definite 3B invcrtible, such that A = BB T (525) 
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A One-dimensional Results 



A.l Gaussian 
A. 1.1 Density 



1 / (x-tf 

p{x) = j= exp [ — 



V2na 2 



2a 2 



A. 1.2 Normalization 



e 2°' 2 as 



I 

J e -(ax 2 +bx+c) dx 
J e c 2 x 2 +c lX +co dx 



V27TCT 2 



exp 



4ac 



-C2 



■ exp 



4a 

c 2 - 4c 2 c 
-4c 2 



A. 1.3 Derivatives 



dp(x) 

dlnp(x) 
dfi 

dp{x) 
da 
d\np(x) 
da 



(x-fi) 



a^ 



p{x)- 



(x-nf 



(x - t lf 



A. 1.4 Completing the Squares 

C2X 2 + c\x + Co = — a(x — b) 2 + w 
, lei 



-a = c 2 



2c 2 



1c? 

w = h c 

4 c 2 



or 



c 2 x + c\x + c 



1 

2T 2 



(x - /i) 2 + d 



-Ci 2 

/" = ^ f = ~ 

2c 2 2c 2 



1 a c ? 

d = c - — 

4c 2 



A. 1.5 Moments 

If the density is expressed by 



p(x) 



I 



\/27rcr 2 



exp 



2(T 2 



(526) 

(527) 
(528) 
(529) 

(530) 
(531) 
(532) 
(533) 



or p(x) = C exp(c 2 x 2 + c\x) (534) 



then the first few basic moments are 
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(x) 

'™2\ 



a 2 + ,i 2 



2c 2 



2c 2 + ( 2c2 ) 



(x 3 ) = 3<r 2 n + yU 3 

(x 4 ) = + 6^ 2 (T 2 + 3a 4 

and the central moments are 

((x-ri) 
((x- M ) 2 ) 
((x- M ) 3 ) 
<(x- M ) 4 ) 



(2c 2 ) 2 



2c 2 



= 

a 2 = 


3a 4 = 



-l 

2c 2 



1 

2c 2 



A kind of pseudo-moments (un-normalized integrals) can easily be derived as 
exp(c 2 x 2 + c\x)x n dx = Z(x n ) 



1 n 


[ c\ 1 


\ exp 


-4c 2 . 


V -C2 



(x"> 



^From the un-ccntralized moments one can derive other entities like 



(535) 



(x 2 > - (x) 2 
(x 3 }-(x 2 }(x> 

(x 4 ) - (x 2 ) 2 



2a 4 + An?u 2 



2 9? 
2ci 

(2c 2 ) 2 

2 

(2c 2 ) 2 



1 -4 



2<; 2 



A.2 One Dimensional Mixture of Gaussians 



A. 2.1 Density and Normalization 



K 



1 (s - A^fc) 5 



(536) 



A. 2. 2 Moments 

An useful fact of MoG, is that 



(537) 



where (•)& denotes average with respect to the k.th component. We can calculate 
the first four moments from the densities 



1 (x - fj, k ) 



p( x ) = E (0fcC ' fc exp [ Cfc2x2 + CfciX 



2 1 



(538) 
(539) 



as 
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( X ) = Efc PkPk 



= Efc Pk 
= J2kPk 

= EfePfc 



2c fc2 



1 l_ | — Cfcl 

i2 \ 2Cfc2 



2c k 



(2c fe2 ) 2 



2c fc2 



If all the gaussians are centered, i.e. p k = for all k, then 
(x) = =0 

( x2 > = EfePfc^fc - EfcPfe 







\ 2 c 2 






(2c M ) 




) - 6 &- 








V 2c <= 2 , 









2c fc2 











( x ) = EfcPfc3c^ = EfePfeS 



n 2 



2c fc2 



^From the un-ccntralized moments one can derive other entities like 



(x 2 )-(a;) 2 = J2k,k>PkPk> 

(x 3 )-(x 2 )(x) = Y,k.k>PkPk> 

(x^)-(xY = E k , k 'PkPk> 

A. 2. 3 Derivatives 



p\ + a t- PkPk>] 

3v 2 kVk + pl-(<y 2 k + pt)Pki] 

4 + 6 M fV 2 + 3ai - (o\ + ^(cr 2 , + p 2 .,)] 



Defining p(s) — J2 k Pk-^siPk, & k ) we get for a parameter 9j of the j.th compo- 
nent 



that is, 



d\np(s) 
89, 



din p(s) 

d Pj 
dlnp(s) 

d P] 
dlnp(s) 



p 3 Af s 3 , o) ) 8 m( P] N s (fij , cr 2 ) ) 

HkPk^siPk^l) 86 j 



HkPk^siPk^l) Pj 

PjK(Pj,<r]) (s-Hj) 



J2 k PkK(p k ,<j 2 k ) 



T,kPkK(p k ,<jl) <Tj 



(s - 



- 1 



(540) 

(541) 
(542) 

(543) 



Note that p k must be constrained to be proper ratios. Defining the ratios by 
pj = e rj I ^2 k e Tk , we obtain 



dlnp(s) _ y , 8 In p(s) 8p t ^crc 



dri 



dpi drj 



Or 



- = Pl{Slj - Pj) 



(544) 
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B Proofs and Details 
B.l Misc Proofs 



B.l.l Proof of Equation 83 



Essentially we need to calculate 
d{X n ) kl d 



dXij dXij 



^ ^ Xk,ui X Ul jU2 X Un _ 1 , 



111 5- ■■ 5^71— 1 

n-1 

^(x'") fcl (x»- 1 -'-), 7 

n-1 



^(X^J^'X"- 1 -'') 



Using the properties of the single entry matrix found in Sec. 9.7.4 the result 
follows easily. 

B.1.2 Details on Eq. [546] 

<9det(X ff AX) = det(X H AX)Tr[(X ff AX)- 1 a(X ff AX)] 

= dct(X H AX)Tr[(X ff AX)- 1 (a(X H )AX + X ff 9(AX))] 
= dct(X if AX)(Tr[(X H AX)- 1 9(X H )AX] 

+Tr[(X H AX)- 1 X H 9(AX)]) 
= det(X ff AX)(Tr[AX(X ff AX)" 1 a(X H )] 

+Tr[(X H AX)- 1 X H Ac)(X)]) 

First, the derivative is found with respect to the real part of X 
adet(X^AX) _ det(xHAX ^Tr[AX(X ff AX)- 1 9(X ff 



Tr[(X H AX)-!X ff A<9(X)] \ 
+ dftX ) 

= dct(X H AX) (AX(X H AX)" 1 + ((X ff AX)^X H A) T ) 
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B PROOFS AND DETAILS 



Through the calculations, ( 92 1 and ( 219 1 were used. In addition, by use of ( 220 1 , 
the derivative is found with respect to the imaginary part of X 



,<9det(X ff AX) 
1 <93X 



idet(X ff AX)( 



/Tr[AX(X /f AX)- 1 a(X H )] 



<93X 

Tr[(X ff AX)- : X ff A<9(X)] \ 



<93X 



det(X ff AX)(AX(X ff AX)- 1 - ((X ff AX)-'X ff A) T ) 



Hence, derivative yields 
<9det(X ff AX) 



1 /<9det(X ff AX) ,<9det(X ff AX) 



ax 



<93?X 



<93X 



= det(X ff AX)((X J/ AX)- 1 X H A) 
and the complex conjugate derivative yields 



<9det(X H AX) 
<9X* 



l^det(X ff AX) ,<9dct(X ff AX) 



dct(X ff AX)AX(X rt AX) 



H 



Notice, for real X, A, the sum of (|228| and ( [229] ) is reduced to (|48l. 
Similar calculations yield 



<9det(XAX H ) 

ax 



1 /<9det(XAX H ) ,<9det(XAX ff ) 



<95RX 



asx 



and 



<9det(XAX H ) 



det(XAX H )(AX ff (XAX H )- 1 ) 



1 /<9dct(XAX ff ) ,<9dct(XAX ff ) 
2 V MX +l dZX 
det(XAX ff )(XAX ff )~ 1 XA 



(545) 



(546) 
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